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ABSTRACT 

Cosmic reionization is expected to be complex, extended and very inhomogeneous. 
Existing constraints at z 6 on the volume- averaged neutral hydrogen fraction, x-^i, 
are highly mo del- dependent and controversial. Constraints at z < 6, suggesting that 
the Universe is highly ionized, are also model-dependent, but more fundamentally are 
invalid in the context of inhomogeneous reionization. As such, it has recently been 
pointed out that there is no conclusive evidence that reionization has completed by 
z ~ 5-6, a fact that has important ramifications on the interpretation of high-redshift 
observations and theoretical models. We present the first direct upper limits on xhi 
at z ~ 5-6 using the simple and robust statistic of the covering fraction of dark pixels 
in the Lya//? forests of high redshift quasars. With a sample of 13 Keck ESI spectra 
we constrain a;Hi^0.2 at 5 < z<5.5, rising to xhi^O.S at z ^ 6.1. We also find 
tentative evidence for a break in the redshift evolution of the dark covering fraction 
at z ~ 5.5. A subsample of two deep spectra provides a more stringent constraint 
of X}ii{z — 6.1) ^0.5 when combined with conservative estimates of cosmic variance. 
This upper limit is comparable to existing results at z ~ 6 but is more robust. The 
results presented here do not rely on assumptions about the quasar continuum, IGM 
density, H II morphology or ionizing background fields, and thus are a good starting 
point for future interpretation of high redshift observations. 

Key words: Galaxies: high-redshift - Cosmology: observations - dark ages, reion- 
ization, first stars - diffuse radiation - early Universe - quasars: absorption lines 



> 1 INTRODUCTION 



The reionization of the Universe by the first generations of 
astrophysical sources is of fundamental importance, offering 
ghmpses into the early Universe and insight into many as- 
trophysical processes. Consequently, much effort has gone 
(and continues to go into) understanding this epoch. Many 
indirect observational probes of reionization have been pro- 
posed. Arguably foremost among these involve the spectra 
of z ~ 6 quasar s, first discovered with the Sloan D igital Sky 
Survey (SDSS: lFan et al.ll200ll . [200l . l2004l2006bl ). 

Quasars have long been useful probes of the IGM across 
a wide range of redshifts through Lya forest studies. The 
moderate redshift {z < 3) Lya forest is largely transparent, 
with occasional transmission gaps arising from strong ab- 
sorption systems, e.g., Lyman Limit Systems (LLSs) and 
Damped Lya Absorbers (DLAs). However, at higher red- 
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shifts (2>5), the Lya forest becomes i ncreasin gly opaque 
and quasar absorp tion spectra saturate ( Becker et al.. .200ll: 
Dior govski et al.l 120011 ) due to the strong absorption cross- 
section of the Lya transition. Even trace amounts of neutral 
hydrogen (with neutral fractions as low as xm ~ 10~^-10~*) 
are then sufficient to render the forest optically thick. Fur- 
thermore, reionization driven by stellar sources progresses 
in a highly inhomogeneous manner, with neutral patches of 
the IGM becoming increasingly rare. Thus, robust and tight 
constraints on reionization from the Lya forest are difficult 
to obtain and require detailed modeling. 

Modeling reionization is a significant challenge, as sim- 
ulations need to be large enough to sta tistically capture 
the biased locations of quasar host halos (|Lidz et al. I l2007l : 
[Alv arez fc Abel 2007,; Mesin ger 2010) as well as the rare 
voids that dominate the transmi ssion at high redshifts (e.g., 
iBecker. Ranch, fc Sargent l2007l ). At the same time, simu- 
lations must have sufficient resolution to resolve the dom- 
inant ionizing population expecte d to reside in atomically 
cooled, M ~ 10*Mq halos (e.g.. iMesinger fc Diikstra,200§ ; 
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IChoudhurv. Ferrara. fc Galleranil 20081). and the small-scal e 



structure in the IGM fe.g.. iLidz. Oh. fc Furlanettd [200^ . 
More approximate techniques have recently been devel- 
oped in order to overcome the challenge s posed by such a 
lajRe r ange of scales, includin R sub-grid (|Kohler fc GnedirJ 
20071: iMcQuinn et al.l l2007bh and semi- numerical (e.g., 



Mesinger fc Furlanettoll2007l ') methods 



Perhaps the most prominent metric to have come out 
of the high-redshift forest studies of reionization is the so- 
called Gunn-Peterson fGP. lGumi fc Peterson|[T965l ) effective 
optical depth, r^p. This quantity is defined somewhat arbi- 
trarily (for lack of an obvious alternative) by the logarithm 
of the mean flux, tqp = — ln{_FGp), and is dependent on the 
complex stru cture of the density, ionization and UV back- 
ground fields l|Fan et al.l2002l : ISongaila fc Cowiell2002l ). Nev- 
ertheless, measurements of Tqp have led to two important 
conclusions: (i) we might be witnessing the final stages of 
reionization, as evidenced by the steep increase with redshift 
in the mean and scatter of tqp at 2: > 6; (ii) the Universe is 
highly ionized at z<&, as evidenced by the lack of such a 
rise and scatter. 

The first conclusion has been the subject of much recent 
controversy. Specifically, empirical extrapolations with alter- 
native (though somewhat arbitrary) models of the density 
field and continuum fitting can reproduce the stee p rise in 
Tqp w ithout h aving to appeal to reionization (Be cker et al.l 
I2OO7I . see also ISongailal [2004 ). Additionally, the observed 
scatter in r^p is consistent with it being driven solely by 
the density field (|Lidz et al.ll2006l '). In fact, it is unlikely that 
the final overlap stages of reionization are accompanied by 
a steep rise in the UV background (or an associ ated fall in 
Tqp) as was initially believed (e.g.. lGnedinll2000l ). since the 
evolution of the ionizing photon mean free path in these end 
stages is likely regulated by L LSs and not reionization itself 
l|Furlanetto fc Mesinge"3l2009l '). 

The second conclusion derived from tqp - that the Uni- 
verse is highly ionized at 2; < 6 - has been far less controver- 
sial. However, similar arguments about the ambiguity of the 
observations at 2 ~ 6 can also be applied at z ~ 5, where 
the standard lore tells us that the Universe has long been 
ionizecO]. Specifically, the photoionization rate, F and mean 
volume- weighted neutral fraction, ih i, are derived fro m rgf, 
assuming a uniform V (e.g., Fig. 7 in iFan et al] |200(j'). How- 
ever, reionization is expected to be highly inhomogeneou^. 



^ Note that Ai: = 1 at these redshifts is a sizable fraction (~ 20%) 
of the Hubble time. Therefore understanding the current con- 
straints on reionization at 2 ~ 5-6 is important not only in cal- 
ibrating reionization models, but also in constraining its im pact 
on subsequent structure formation (e.K., |Busha et al.[l201Cl ). and 
in interpreting virtually all high redshift observations. 
^ Reionization driven by stellar sources is "inside-out" 
on large-scales: ionized bubbles grow around the first, 
highly clustered sources, with their eventual coalescence 
( "overlap") signaling the completi on of reionization (e.j 



Furlanetto. Zaldarriaga. fc Hernguis^ |2004| : iMcQuinn et al.l 
2007bl: iTrac fc Cer |2007|: IZahn et al.ll201(]| : see the recent review 
by iTrac'^T'Cnedin |2009| ). These reionization scenarios result 



in a very inhomogenous ionization field. If on the other hand, 
quasars (or other ionizing sources with a hard spectrum) drive 
the bulk of reionization, the ionization fields would be more 
uniform. However, the abundance of quasars seems too low for 
them to contribute a significant fraction of ionizing photons at 



with pre-overlap neutral regions having a very weak ioniz- 
ing background (and hence low F). Therefore, the a priori 
assumption of a uniform F is tantamount to assuming that 
reionization has already completed; using the derived F or 
iHi as evidence of this is circular logi(0. The Lyman forests 
can be made opaque by both (i) regions with trace amounts 
of neutral hydrogen inside the ionized IGM (shi ~ 10~^- 
f0~*), and (ii) the neutral patches pre-overlap (xhi ~0.1-1). 
Observationally, we cannot distinguish between these. The- 
oretically, uncertainties in the relevant parameters are large 
enough to allow a sizable contribution from (ii), i.e. aSni ^ 
tens of percent. 

Given these complications, how best can existing quasar 
spectra be used to constrain reionization (i.e. shi), in the 
least model-dependent fashion? The strength of the Lya and 
Ly/3 transitions is such that even trace amounts of neutral 
hydrogen will suppress the observed flux well below the de- 
tection limits of current instruments, thus any remaining 
neutral patches along the lines-of-sight to distant quasars 
will generate "dark" pixels in the spectra of those quasars. 
Mesinger (2010) argued that the simple metric of the cover- 
ing fraction of these dark spectral patches provides the least 
model-dependent constraint on the neutral fraction. This 
constraint is an upper limit since it does not differentiate 
between the above-mentioned (i) and (ii). At z>5, the con- 
tributions to absorption from the thickening Lyman forests 
and from the evaporating residual H I patches will be ex- 
ceedingly difficult to disentangle without detailed modeling 
of the IGM. Even with complete model confidence, the cur- 
rent observational sample is insufficient to constrain xm to 
below the percent level (,Mesingcr 201Q ). 

In this work, we ignore these complications and instead 
present a simple upper limit on the neutral fraction at z ^ 
5-6 using the covering fraction of dark pixels. While this 
limit is highly conservative since it does not "model-out" 
the substantial absorption from residual H I in the ionized 
IGM, we stress that this is the only observational constraint 
derived from existing quasar spectra free of a priori, model- 
dependent assumptions. As such it is a natural first step 
towards interpreting the observations. 

This work is organized as follows. In ij^] we describe our 
observational sample of 2 ~ 6 quasar spectra, while in ij3]we 
discuss our methodology for deriving the covering fraction 
of dark pixels from this sample. In ij?) we show example 
spectra. In ^we present the corresponding upper limits on 
xhi, which are the main result of this work. In iJO] we make 
simple estimates of the sightline-to-sightline variance dur- 
ing reionization, and in [J7]we discuss future work. Finally, 
we summarize our results in iJS] We quote all quantities in 



2 > 5 l lFan et al.l I2OOII : Ijiang et al.l |200S| : IWillott et afl |201C|) : 

furthermore, even if faint quasars missed by existing surveys 
were present in sufficient abundance, t hey would overproduce the 
unres olved soft X-ray background l lDiikstra. Haiman. fc Loebl 
12004) . Therefore, even in extrem e scenarios involving subst antial 
pre-ionization by quasars (e.g., IVolonteri fc GnedinI [20091) , the 
late stages of reionization are expected to be inhomogeneous. 
^ Despite the sizable spatial fluctuations in F post-reionization, 
averaged statistics such as r^p are still dominated by the spatial 
fluctuations in the density field. Therefore the assumption of a 
uniform F is valid post-reionization, since it only biases the in- 
ferred value of F by a few percent ( Mesinger fc Furlanetto 20091) . 
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Table 1. Keck ESI spectra 



Object 




z 




tcxp (hr) ''"GP,li 


tn 


Rcf. 


J0927+2001 


5 


77 


19 


88 





33 


3 


2 


1,8 


J0836+0054 


5 


81 


18 


74 





33 


4 





2 


J0840+5624 


5 


84 


19 


76 





33 


3 


4 


1,8 


J1335+3533 


5 


90 


20 


10 


Q 


33 


3 


1 


1,8 


J1411+1217 


5 


90 


19 


63 


1 


00 


3 





3,8 


J0841+2905 


5 


98 


19 


84 





33 


2 


9 


4 


J1306+0356 


6 


02 


19 


47 





25 


3 


7 


2,8 


J1137+3549 


6 


03 


19 


54 





67 


3 


2 


1,8 


J0353+0104 


6 


05 


20 


54 


1 


00 


2 


9 


5 


J0842+1218 


6 


08 


19 


64 





67 


3 


4 


6 


J1623+3112 


6 


25 


20 


09 


1 


00 


3 


6 


3,8 


J1030+0524 


6 


31 


20 


05 


10 


32 


4 


6 


2,8,9 


J1148+5251 


6 


42 


20 


12 


11 


00 


5 


1 


7,8,9 



Refer ences: 1) iFan et al]l2006tl: 2) iFan et al]|200ll: 3) iFan et al l 
2004 4')[G oto 260( 
7) iFan et al...2003 : 



|2004 4) [G oto 200^; 5) I Jiang e t aL 200'^ ; 6) Ide Rosa et al|l201ll : 
) iFan et al.,.200&: 9) IWhite et al.ll2003l . 



comoving units using a standard ACDM cosmology with pa- 
rameters (^A, ^M, fit, n, as, Ho) = (0.72, 0.28, 0.046, 0.96, 
0.82, 70 km s~^ Mpc~ ^), consistent with re c ent results from 
the WMAP satellite (|Dunklev et all l2009l : iKomatsu et all 
I2OO9I '). 
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Figure 1. hya/p redshift paths spanned by the ESI spectra. For 
each quasar the redshift is indicated by a star, the Ly« forest by 
the upper (blue) line, and the Ly/3 forest by the lower (red) line. 



2 THE DATA 

Our sample consists of 13 quasars at z ~ 6 with spec- 
troscopy from the E chellette Spectrograph and Imager (ESI; 
ISheinis et all 120021 ') on the Keck II telescope (Table [TJ. We 
use only ESI spectra in order to have uniform wavelength 
coverage and resolution, and also for the high signal-to-noise 
ratio (S/N) even in relatively sho rt exposures. Details of the 
reduction procedures are given i nlWhite et al.l ( 2003). Most 
of these spectra were included in lFan et al. ( 20061 ) : as in that 
work we exclude quasars with broad absorption line (BAL) 
features. The quasa r redshifts ha v e been updated to match 
the values given in ICarilU et al.1 (|2010l ). taking advantage 
of more accurate CO and Mg II redshifts when available. 
The spectra have widely varying exposure times (from 15 
minutes to > 10 hours) and span 1.8 optical magnitudes in 
flux, and thus have widely differing S/N levels in the Lyman 
forest regions. 

We quantify the dynamic range of each spectrum in 
terms of the average limiting optical depth in the Lya and 
LyP forests. To derive these quantities, we fit t he power-law 
contin uum model for the rest-frame UV from iTelfer et al.l 
l|2002l) to the unabsorbed continuum redward of Lya. The 
average 2a limiting optical depths in the Lya forest are then 
calculated as Tqp — — (ln(2cr//cont)) and shown in Ta- 
ble [T] For the Ly/3 forest the continuum flux, and hence 
dynamic range, is affected by absorption from the under- 
lying z ~4-5 Lya forest, as the photons redshift into Lya 
resonance at 1 -|- = (A^/Aq)(1 -|- Zjs). We account for this 
eff ect using the red shift evolution of the mean opacity given 
m iFan et al.1 (|2006l ). On the other hand, the optical depth 
is reduced by the weaker oscillator strength (/) and shorter 
wavelength of Ly/3 (tgp oc /A), which increases the dynamic 
range in the Ly/3 forest relative to Lya. Note that the con- 



tinuum fitting and the derived optical depths are strictly 
used to characterize the quality of the spectra, and do enter 
into the classification of "dark" pixels, which is the main 
result of this work (see 

We bin the spectra into pixels of width 3.3 Mpc. The 
ESI spectra have a wavelength resolution of i? ~ 3000 — 6000 
(IFan et al. 2006 ) , though they are sampled at a higher reso- 
lution during extraction. Binning to 3.3 Mpc results in pixels 
with R ~ 800 calculated from the mean of ~ 30 unbinned 
pixels. This increases the dynamic range by ~ 20 — 50% (as 
determined from the optical depth) 



* Our choice of pixel size is fairly arbitrary, and could be con- 
sidered one of the few "model" uncertainties in this work. This 
scale is a factor of a few times the Jean's length in the mean- 
density, ionized IGM. Furthermore, a simple estimate of the pho- 
toevaporation timescale of me an-density, -^3 Mpc regions exposed 
to the observed mean F (e.g. jBolton fc Haehneiil2007d) , yields 
Az ~ 0.8. This suggests that even smaller H I patches could 
persevere for a non-trivial time. Neglecting such sub-pixel H I 
structure likely makes our results again more conservative, since 
the associated strong opacity in the core of the Lyman line, as 
well as a potentially strong contribution from the damping wings 
fe.g.. lMiralda-Escud6ll998i i. would likely cause the entire 3.3 Mpc 
"pixel" to be dark, resulting in an overestimate of the covering 
fraction. More generally, by absorbing flux in neighboring pix- 
els, damping wing absorption from potential H I patches would 
make our upper limits more conservative. Unfortunately, estimat- 
ing the impact of damping wings is highly model-dependent, re- 
quiring the distribution of H I impact parameters and partial 
ionization levels. In either case, we plan on combining our sky- 
noise-limited ESI sp ectra with read-noi se-limited high-resolution 
HIRES spectra (e.g.. iBecker et al.|[2007l l in an attempt to resolve 
sub-structure. Note that we are not discussing collapsed struc- 
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Though the ESI spectra have been degraded in resolu- 
tion in order to increase the dynamic range, the higher native 
resolution of the spectra confines the strong night sky emis- 
sion lines to relatively small scales. We take advantage of 
this property to select bins between the brightest sky lines. 
First, the average noise is calculated in a 3.3 Mpc window 
centered on each pixel in the unbinned spectrum. Minima in 
the average noise are identified, and bins are placed at these 
minima. This effectively identifies locations where binned 
pixels can "fit" between bright sky lines. Finally, the spaces 
between these bins are iteratively filled with additional bins 
located at the minima within the spaces until no more bins 
can fit. We then remove from our analysis any pixels with 
too low of a threshold in optical depth, defined as a 2a limit 
of Ta = 2.5. 

This process results in relatively uniform noise across 
the full wavelength range, without resorting to smoothing 
techniqu es that could amp lify correlations between the pix- 
els (e.g-. IWhite et al]|2003l ). In addition, most of the binned 
pixels are non-adjacent, further reducing the susceptibility 
of our calculations to correlated noise introduced during ex- 
traction of the 2D spectra, and also to any physical corre- 
lations within the IGM on scales corresponding to the bin 
size. In the following sections we use "pixel" to refer to a 3.3 
Mpc binned pixel derived from this process. 

Finally, we limit the wavelength ranges of the Lyman 
forests used in our analysis in order to avoid biases and 
contamination. The red edge of both the Lya and Ly/3 
forests is taken to be 40 Mpc from the quasar redshift 
(A2; = 0.1). This is roughly the minimum distance required 
to reduce bias from the quasar's surrounding H II morphol- 
ogy due to the local oyerabunda.nce o f clustered galaxies 
l|Alvarez fc Abell l2007l : iMesingeJ I2OI0I ). The blue edge of 
the Lya forest is at 1 -I- z'^i-^ = (1 -I- 2;em)(1040A/AQ), cho- 
sen to avoid wave lengths affected by emission from Ly/? and 
O VI (|Fan et al.l [2OO6 I , and contamination from the Ly/? 
forest. The blue edge of the Ly/3 forest is at 1 -I- z^^^ — 
(1 -I- Zem)(970A/A^) in order to avoid contamination from 
Ly7 emission and the Ly7 forest. We show the span of the 
Lya and Ly/3 forests in Fig. [1] for all spectra used in our 
analysis. 



3 DEFINING "DARK" PIXELS 

We now quantitatively define what we mean by a "dark" 
pixel. We apply two criteria, as defined below, to the Lya 
and Ly/J forests from our spectra. In addition, we note that 
dark pixels arising from a pre-overlap neutral region must 
be dark m both Lya and Ly/?. Therefore, where we have 
overlapping (in redshift, see Fig.[l]) spectral coverage of Lya 
and Ly/3 we also present the combined dark fractions. 

We stress that neither method used to calculate the cov- 
ering fraction of dark pixels directly relies on continuum fit- 
ting, and thus we are relatively insensitive to the significant 
uncertainties associated with this process. 



tures (e.g.. Illiev. Shapiro, fc RaEall2005l ). which are more appro- 
priately labeled as LLSs or DLAs contributing to the absorption 
inside the ionized IGM component. 



3.1 Dark pixel fraction using a flux threshold 

First, we define dark pixels using only the flux and noise 
properties of the spectra: a pixel is dark if it has fiux 
fi < noi, where fi is the measured flux and ai is the mea- 
sured noise of each pixel. Simply put, a pixel is dark if it has 
a flux consistent with zero within the given (assumed Gaus- 
sian) noise, and is not dark if it has {S/N)i > n. We have 
adopted a threshold of n = 2. For this threshold, ~ 2.3% of 
dark pixels will be scattered above the threshold and mis- 
classified; thus we augment our dark pixel counts and un- 
certainties by this amount. Additionally, pixels with fiuxes 
just above the threshold can be scattered below it, making 
our upper limits slightly more conservative (modeling this 
minor effect would require knowledge of the underlying pixel 
fiux distribution). 

Our choice of a 2a threshold is somewhat arbitrary. We 
compared results obtained with a 3cr threshold and found 
they are consistent - though understandably slightly weaker 
- with those presented below. 

3.2 Dark pixel fraction using negative pixels 

If the noise in the extracted ESI spectrum is symmetric 
about zero, pixels that intrinsically have zero fiux will have 
equal probabilities to be scattered into positive and nega- 
tive values. Therefore, we can also estimate the number of 
dark pixels by simply counting pixels with negative flux and 
multiplying by two. Again, some pixels with small positive 
fiuxes will be scattered to negative fluxes and thus incor- 
rectly counted as dark, but this effect should be even smaller 
than when using the flux threshold method. 

When we present the combined dark fraction from 
redshift-overlapping Lya and Ly/3 forests, we multiply the 
number of instances of negative flux in both Lya and Ly/3 by 
a factor of four. Again, if the noise is symmetric about zero, 
then a pixel with intrinsic zero-ffux will have a 1/4 chance 
of having observed negative ffuxes in both Lya and Ly/3. 

Counting dark pixels using negative values relies on 
the assumption that the noise is symmetrically distributed 
about zero for pixels with no intrinsic flux. Systematic effects 
in the reduction process, such as incorrect sky subtraction, 
could bias these values and hence our results. We do not see 
any evidence for this type of bias from examination of the 
pixel flux distributionQ Additionally, the deeper spectra are 
less susceptible to this type of bias, as they are created by 
combining many individual exposures, averaging over such 
effects. 



4 LYa//3 FOREST EXAMPLES 

Figure [2] shows example Lya forests for two objects. The 
figure demonstrates the placement of bins away from bright 
sky lines and the relatively uniform noise in the binned pix- 
els used for the final calculation of the dark pixel fraction. 
The forest for J0841-(-2905 is considerably darker than that 

^ For example, examination of the flux distributions in the forest 
regions normalized by the standard error shows that they arc 
roughly Gaussian with fi = and cr = 1, ignoring the tail to 
positive flux values arising from flux leaks within the forest. 
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5.6 

redshift [Lya] 

Figure 2. Lyo forests of J0841+2905 {zq = 5.98, top panel) and J1306+0356 (zq = 6.02, bottom panel). Pixel flux measurements are 
indicated with bars, and shown in gray when flux is detected (according to the 2a threshold), and black for dark pixels. Dark pixels with 
fiux ^ are indicated with small squares at negative fluxes. The smooth red line shows the 2(t noise level and defines the dark threshold. 
The orange line shows the unbinned noise, scaled to match the la binned noise level. This highlights the placement of the bins between 
the stronger sky lines. Finally, the dashed magenta line shows the fitted continuum scaled by e~^, with t = 2.5. Pixels are rejected when 
the binned noise exceeds this threshold, as in the case with J084f+2905 for three pixels at z > 5.8 (shown as hatched bars). 



5.45 5.50 5.55 5.60 5.65 5.70 



10 
8 

" 6 

i/i 4 

ai 

t 1 
S 2 
< 3 
4 
5 
6 



m 







Ly/3 



"B- 



5.45 5.50 



5.55 5.60 
redshift [Lya//3] 



5.65 5.70 



Figure 3. Comparison of the Lya and Ly/3 forests for 
J0836+0054 aligned in absorption redshift (note that the flux 
scale for the Ly/3 pixels is inverted for easier comparison). The 2a 
noise level is shown with a red line. The excess of dark pixels in 
the forest is likely due to contamination from the overlapping 
Lya forest at 1 + = {Xp / Xa){l + z/j). The pixels are aligned 
in redshift. Comparing the two forests shows that only two pixels 
are dark in both. 



of J 1306+0356 in the same redshift range. This is hkely a 
consequence of the lower dynamic range of the J0841+2905 
spectrum, but there could also be contributions from line- 
of-sight variance. This highhghts the need for both many 
lines-of-sight and high quality spectra when deriving con- 
straints from the Lya forest. 

The Ly/3 forest yields a stronger constraint on the frac- 
tion of dark pixels, as the relative weakness of the Ly/3 tran- 
sition compared to Lya results in greater dynamic range. 
The specific ratio of t"" jr^ varies from pixel-to- pixel de- 
pend i ng on the sub-pi xel structure of the IGM (e.g.. lSongailal 
I2OO4I : [Fan et al.ll2006l ). and the absorption from the under- 
lying Lya forest (see Since we are treating the dark 



fraction as a direct upper limit on Shi, we do not introduce 
any model assumptions about either the dumpiness of the 
IGM at 2 ~ 5-6 or the nature of the lower-redshift Lya 
forest into our results. In future work, we will attempt to 
resolve and model-out the underlying Lya forest from the 
Ly/3 region, yielding tighter constraints, albeit at the cost of 
introducing model uncertainties (see 

In Fig. |3] we show a comparison of the Lya and Ly/3 
forests for J0836-I-0054, overlayed in absorption redshift. 
The 2a noise level is shown with a red line, and dark pixels 
identified with this flux threshold (i3]T| are shown in black. 
As mentioned above, the contamination of the Ly/3 by the 
lower-redshift Lya forest can result in dark pixels even when 
the corresponding pixels in Lya are not dark. Pixels corre- 
sponding to pre-overlap, highly neutral regions must be dark 
in both Lya and Ly/3. In this figure there are two such pix- 
els. The intersection of dark pixels in both forests provides 
the strongest possible constraint available from this method; 
however, due to the limited dynamic range the presence of 
dark pixels in both forests is a necessary but not sufficient 
indication of neutral IGM. 

Finally, it is important to remember that each of our 
pixels is a non-trivial ~ 3.3 Mpc in size. As already dis- 
cussed, a single dark pixel is sufficiently large to be a poten- 
tial tracer of pre-overlap HI, especially when weighted by the 
LOS impact parameter. In fact, extended dark patches need 
not be a good tracer of reionization in high-reds hift spec- 
tra, e specially when the spectra begin to saturate (jMesingeil 
I2OIOI ). 



5 RESULTS 

The main results of this work are presented in Fig. 3] The 
covering fractions of dark pixels computed using the flux 
threshold (i3TT]) and negative pixel ( ^3.2|l statistics are di- 
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Figure 4. Upper limits on the neutral hydrogen fraction derived from Keck ESI spectra of ^ ~ 6 quasars using Lya and Ly/3 forest 
absorption. The constraints are obtained from the covering fraction of dark pixels in redshift bins of width Az = 0.15 (Lyo) and Az = 0.2 
(Ly/3). The upper left panel shows the number of pixels in each redshift bin, and the upper right panel shows the number of independent 
lines-of-sight contributing to each bin. The lower left panel shows the constraints derived by counting dark pixels using a 2a flux threshold 
( i|3.1ll : the lower right panel shows the constraints derived from pixels with negative fluxes ( i|3.2ll . The number of pixels and lines-of-sight 
are the same for both methods. Uncertainties are calculated from jackknife statistics (see text), as we are only concerned with an upper 
limit the uncertainty is only shown on the upper bound. 
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redshift 
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redshift 



6.2 



Figure 5. Same as Figure[4] but using only the subsample of our two deepest z > 6.2 spectra (> 10 hour total exposure times). With only 
two spectra, these constraints are more susceptible to line-of-sight variance. Instead of the Jackknife errors in Fig. (4] we plot estimates of 
the uncertainty due to the cosmic variance of reionization at zhi ~ 0.1 (see discussion in These error bars, although conservatively 
large, are admittedly model-dependent. 



rectly converted into a upper limil[f| on the neutral hydrogen 
fraction and shown in the bottom left and rigfit panels, re- 
spectively. The dark fractions calculated from the Lya, Ly/3, 
and combined Lya -I- Ly/3 forests are shown in blue, red, 
and magenta, respectively. Pixels are binned in redshift bins 
of A« ~ 0.15 (0.2) for the Lyo (Ly/3) foreslQ. The upper 

^ As mentioned previously, in addition to any possible pre-overlap 
H I, the dark pixel fraction also includes contributions from H I 
embedded in the ionized IGM with an optical depth greater than 
the available dynamic range; hence it is only an upper limit on 
Shi- Additionally, partially ionized neutral IGM (see discussion 
in Q and unresolved ionization substructure (on scales smaller 
than our pixels) make our upper limits even more conservative. 
^ Wide redshift bins decrease the sample variance, but also aver- 
age over any IGM evolution. Our fiducial bin choice (although 
fairly arbitrary) is narrower than the likely reionization time- 
scales: Hubble time, photoevaporation time of pixel-sized re- 
gi ons, and the growth of the halo mass function (see, e.g.. Fig. 9 
inlMesinger, Johnson, fc Haimaiill2006l and Fig. 1 in iLidz et al.l 
I2OO7II . Therefore, our bins are likely narrow enough to resolve 



left panel shows the number of pixels in each redshift bin, 
and the upper right panel shows the number of independent 
lines-of-sight contributing to each bin. 

The fraction of dark pixels in each bin and the asso- 
ciated uncertainty on this fraction are obtained using the 
jackknife method. We repeat the calculation 13 times, each 
time leaving out one of the spectra. We then use jackknife 
statistics to derive the mean and variance from the fractions 
calculated during each resampling. Since the dark fraction 
is an upper limit on xm, the derived variance is shown only 
in the upper limits of Figure [J] Jackknife errors may or may 
not fully capture the variance inherent to the reionization 



cosmic evolution, and wide enough to have decent sampling of 
the (possible) reionization morphology in the xhi^O.I regime 
(see 

* Note that this is not a strict application of the jackknife 
method, as leaving out one spectrum at a time results in a vary- 
ing number of pixels in each resampling. However, this variation 
is small. 
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process; in any case, we estimate such cosmic variance to 
be small for our full sample (see !j6]). Rather, the sample 
variance is most likely dominated by the disparate dynamic 
ranges of the spectra. 

Using the flux threshold (negative pixel) method, we 
find upper limits of xhi<0.4 (0.2) at 5.0<2<5.5 (fairly 
independent of redshift), increasing to xm < 0.7 (0.7) at 
2 ~ 6. As expected, the most powerful constraint for the 
flux threshold method is obtained from the combined Lya 
and Ly/3 forests, by requiring that overlapping pixels in both 
forests are darfl The benefit of using the overlapping forests 
is not as evident for the negative pixel method, since it re- 
lies on events rarer by a factor of two (see H3.2p , and so has 
intrinsically larger sample variance than using just a sin- 
gle forest. Including the 1-cr jackknife errors degrades the 
z > 5.5 constraints by ~ 10% (~ 20%) for the flux threshold 
(negative pixel) method. 

The two methods for estimating the dark covering 
fraction yield comparable results, though the negative flux 
method generally results in lower dark fractions. This could 
be due to many pixels with small positive flux levels con- 
tributing to the flux threshold method counts ( jjS.lf) . By 
pushing the threshold to flux < many such pixels are elim- 
inated, at the cost of higher sample variance. 

There is a notable change in the evolution of the dark 
fraction at z ~ 5.5. This feature appears to be genuine and 
is due to the fact that dark gaps in several spectra terminate 
at this redshift. Although it might be tempting to suggest 
this feature could be due to reionization, such an interpreta- 
tion is premature. There is no obvious, corresponding feature 
in the redshift evolution of r^p (|Songaila 2004; Fan et al., 
l2006l : iBecker et al]|2007l ) . Nevertheless, given the large scat- 
ter in the tqp sample, and the difficulties in using tqp as a 
reionization tracer, it will be worthwhile to investigate this 
feature with larger sample sizes and modeling. 

In Fig. O we present results using only our two deep- 
est spectra (>10 hour exposure times). These have an ef- 
fective dynamic range of tqp jj^^ ~ 5 and r^p ~ 7-20 
(see Table With just the two deepest spectra, we ob- 
tain mean dark fractions of xm ^ 0.3 (0.2) at 2: < 5.5 and 
xhi ^ 0.85 (0.5) at 2 ~ 6, with the flux threshold (negative 
pixel) criteria. 

Although these mean values are generally lower than 
those from our full sample (Fig. U), they are also more sus- 
ceptible to cosmic variance. Having two LOSs increases the 
fractional uncertainty in the sample mean by a factor of 
~ 2.5 compared to our full sample of 13 LOSs, since the 
uncertainty in the sample mean scales with the number of 
samples as 1/ V A^los (see |6}. This susceptibility to cosmic 
variance depends on Shi, increasing as the neutral patches 
become rarer. Taking the speciflc case of Shi ~ 0.1 from Fig. 
[6l we see that the uncertainty in the mean dark fraction (due 
to reionization) from two LOS segments of width Az = 0.2 
is comparable to the mean. Reducing this uncertainty re- 



^ Note that the bins used in the analysis of the Lya, Ly/3, and 
Lya -I- Ly/3 forests as shown in Fig.|4]are offset in redshift. There- 
fore the sample variance between the bins can result in a Lyo -|- 
Ly/3 dark fraction which is slightly higher than one obtained just 
from the Ly/3, as is the case for the last set of points in the lower 
left panel. 



quires more deep spectra in order to take full advantage of 
the increased dynamic range. 

Even with this model-dependent, though conservative, 
choice of cosmic variance error bars, our deepest two spec- 
tra yield xm ^ 0.3 at 5 < 2 < 5.5, comparable with those ob- 
tained from the full sample. Furthermore, using the negative 
pixel criteria in the right panel of Fig. [S] we are able to place 
a robust upper limit at z ~ 6.1 of Shi ^0.5. 

Several model-dependent constraints on reionization 
at z ~ 6 have been derived from other astrophysi- 
cal probes such as: (1) the size of the prox i mity zone 
around quasars (IWvithe. Loeb. fc Carillil l2005l: iFan et al 



20061 : ICariUi et aLll201Gl. but se e iMesinger. Haiman fc Cen 
20o4 iBolton fc Haehneld l2007l : iMaseUi et al.l l2007h : (2) a 



claimed detection of damping wing absorptio n from neu- 
tral IGM in quasar spectra (|Me singcr & Hai manI |2004| . 
l2007l . but see iMesinger fc Furlanettd 2008); (3) the non- 
detection of intergalactic dampi ng wing ab s orptio n in 
a g amma ray bur s t spe ctrum (|Totani et al.l |2006| , but 
see iMcOuinn et al.l I2OO8I ): and (4) the number density 



and clustering of Lyg emitters (iMalhotra fc Rhoadsl 2004 
Hai man fc Cenll2005l: iFurlanetto. Zaldarriaga. fc HernquistI 



200 d: iKashikawa et all I2OO6I: iMcQuinn et al.ll2007l. but see 



Dii kstra. Wvithe. fc Haimar] 20071 : Mesinger fc Furlanettol 
i2008bl : llliev et al 1 12008| . Diikstra. Mesinger, fc Wyithe, in 
preparation). All of these constraints are controversial, with 
considerable uncertainties. Our limit of xhi{z — 6.1) < 0.5 
is comparable to or better than these existing upper limits, 
and is much more robust (the model-dependence comes in 
the form of the cosmic variance error bar, and our choice 
here is conservative). 

At lower redshifts {z ~ 5-6) , previous claims of a highly- 
ionized IGM (Shi ~ lO'^'-lO"'*; e.g.. Fig. 7 in iFan et all 
2006) are derived under the a priori assumption of a uni- 
form background. Given that reionization by stellar sources 
is highly inhomogenous, this carries the implicit assump- 
tion that reionization is over, and therefore the standard 
Tqp — > Shi conversion cannot be used to place direct con- 
straints on reionization, even with complete confidence in 
the density distribution of the IGM. The upper limits we 
present in Fig. |4] are the first direct, model-independent con- 
straints on Shi a,t these redshifts. 



6 COSMIC VARIANCE DURING 
REIONIZATION 

The jackknife error bars in Fig. 3] show the sightline-to- 
sightline standard deviation in the dark pixel fraction mea- 
sured from our observational sample. The main source of 
this variation is the disparate dynamic range of our spectra. 
However, the covering fraction of dark pixels itself can in- 
trinsically vary from sightline-to-sightline. Again, this dark 
fraction has a contribution from the ionized IGM and poten- 
tially the pre-overlap neutral IGM. Since this paper places 
upper limits on Shi, we concern ourselves with the cosmic 
variance of the latter, i.e., reionization. 

During the final stages of reionization, the pre-overlap 
neutral regions become increasingly rare. If the mean spac- 
ing of the neutral islands is comparable to or greater than 
the length of our LOS segments, then our sample is prone to 
cosmic variance, and we would need many LOSs to have ro- 
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Figure 6. Fractional uncertainty, cr/fi, of the mean covering frac- 
tion of neutral IGM obtained from A^los samples (see text for de- 
tails). Curves correspond to xm = 0.25, 0.1, 0.014, 0.0013, bottom 
to top. The curves indicate how many Az = 0.2 LOS segments 
are required for their sample average to be accurate within the 
given fractional uncertainty. The shaded region corresponds to the 
sample used in this work. The dashed vertical line roughly corre- 
sponds to the current available sample of 2: > 5 quasars (though 
most currently do not have deep spectra and their forests do not 
all overlap). All curves are computed at 2; = 5 and assume LOS 
segments of length Az = 0.2. 

bust constraints on xm- This intrinsic cosmic variance might 
not be included in empirically derived jackknife error bars. 
To lessen the cosmic variance, one could extend the size of 
the LOS segments (i.e. redshift bins), but then one is aver- 
aging over redshift evolution. 

Here we try to get a rough, quantitative estimate of the 
number of LOSs required for an accurate sampling of pre- 
overlap H I. We note that this is a highly model-dependent 
estimate and so a detailed treatment is not appropriate for 
this work. We begin by noting that the covering fraction 
of neutral IGM during reionization, computed from an en- 
semble of LOS segments, can be thought of as sampling 
some fundamental distribution with mean and stan- 
dard deviation, ctilos- According to the central limit the- 
orem, as the sample size A^los increases, the sample aver- 
age of the covering fraction of neutral IGM approaches a 
normal distribution with mean ^ and standard deviation, 
<y = (TiLOs/V A'los- 

In Fig. |6l we plot the fractional uncertainty, a/fi, of the 
mean covering fraction of neutral IGM obtained from A^los 
samples. The standard deviation of the single sample (i.e. 
from 1 LOS) distribution of the H I covering fraction, ctilos, 

Note that since the "neutral" IGM during reionization could 
be partially ionized, the mean covering fraction, fi, does not have 
to equal xhi, but can be larger. Because of the strength of the 
Lyman series transitions, the covering fraction statistic does not 
distinguish between partially ionized H I regions, 3;hi ^ 0.01, and 
those which are fully neutral, xjii ~ 1. However, the mean neutral 
fraction, xhi, does care if regions are partially or fully ionized. 
Note that here we are only referring to the "neutral" IGM, not 
the residual H I in the ionized IGM, whic h sh ould be at the level of 
< 10"'' (e.g.,|^lton & Hachnelt 2007b|). For example, in the 
fiducial model of iMesi ngerl 1(20101), for xhi = 0.1 we obtain fj, = 
0.18. Therefore, we differentiate between the covering fraction of 
neutral IGM and xm. We caution the reader that this difference 
is highly model dependent. 



is taken from the fiducial models of iMesingej |20l3E3- AH 
curves are computed at z — 5 and assume LOS segments of 
length Az = 0.2, roughly corresponding to the bin size we 
use in this work. Roughly speaking, the curves indicate how 
many LOS segments are required such that their average 
covering fraction is close to the true value within the given 
fractional uncertainty. 

The shaded region in Fig. [6] corresponds to the sam- 
ple used in this work (see the top, right panel of Fig. |4]). 
This means that at 2 < 5.8, where we have a sample size 
of A'los ~ 10, we can be reasonably certain that our sam- 
pling of pre-overlap H I is accurate to better than ~50% of 
the true value, if the true value is ^hi^O.I (i.e., a signifi- 
cantly neutral Universe). Of course, if the Universe is more 
ionized, the neutral IGM patches are rarer, and thus more 
LOSs are needed to obtain the same fractional uncertainty 
in their sample mean. For example, obtaining comparable 
accuracy for xm ~ 0.01 requires ~ 100 LOSs, which is just 
out of reach of the current published set of 2 ~ 5 spectra, 
given that their forests do not all overlarF^. However, when 
the neutral patches are that rare, the dark covering fraction 
statistic is dominated by the ionized I GM, and is a very poor 
probe of reionization (|Mesingeij|2010l ). 

Constraints derived from our two deepest spectra (Fig. 
[5| are more susceptible to cosmic variance than those de- 
rived from our full sample. The fractional uncertainty in the 
mean value of the covering fraction of neutral IGM is a fac- 
tor of 2.5 higher for two of our spectral segments than for 
13. From Fig. [SI we see that for two spectra at xm = 0.1, the 
standard deviation of the sample average is a factor ~ 0.8 
times the mean value of /i = 0.18, i.e. a = 0.8 x 0.18 — 0.14. 
For illustrative purposes, we use these estimates of the cos- 
mic variance uncertainty for our dark fraction points in Fig. 
[5] Note that these error bars are model-dependent theoret- 
ical estimates, unlike the empirical jackknife errors shown 
for our full sample in Fig. |4l where cosmic variance was less 
problematic. 

The error bars in Fig. [5] are appropriate for the covering 
fraction of neutral IGM, which can be substantially higher 
than xhi, due to partial ionizations. In keeping with the 
conservative spirit of this work, we do not convert from the 
covering fraction to xui even in the error bars, making them 
overestimates of the uncertainty on xhi- 



7 FUTURE WORK 

The simple upper limit on xm obtained from this study is 
limited by contamination from absorption within the ionized 



These reionization morphologies are created with 
the publicly available, semi-numerical simulation DexM: 
|http : //www ■ astr o ■ pr incet on . e du/ - mes inger / S im . html Simu- 
lation boxes are 2 Gpc on a side with 3 Mpc cells. In the fiducial 
models, suites of reionization morphologies are created at a fixed 
redshift by varying the ionization efficiencies of atomically-cooled 
(Mniin = 2 X 10^ Mq) sources, with an assumed attenuation 
length inside the ionized IGM of 50 Mpc. 

Of course constraining xm to percent level precision with the 
Lyman forest is difficult to imagine, and so this discussion is 
mostly academic. 
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IGM. This contamination can be addressed on two fronts: 
improved observations and theoretical modehng. 

The observations can be improved by increasing the dy- 
namic range of the pixels in order to eliminate pixels with in- 
termediate absorption levels. Figure [5] shows that very deep 
spectra can strengthen the xui constraint by a factor of two; 
however, obtaining deep spectra is highly costly in telescope 
time, and the payoff in dynamic range with exposure time 
rises only as 5t oc 0.51n(t). In addition, the Ly/3 forest pro- 
vides a stronger constraint than Lya, but has a much shorter 
redshift path, and thus the amount of overlapping Ly/3 for- 
est between quasars at a range of redshifts grows slowly with 
the number of objects (note there is only one redshift bin for 
Lya+Ly/3 in Fig. [5]). An intermediate approach of obtainin- 
ing a few hours of integration time on many objects seems 
better: balancing the needs for depth to detect intermedi- 
ate absorption and for many independent lines-of-sight to 
reduce cosmic variance. This program can be accomplished 
with the sample of ~ 20 bright (zab < 20.5) quasars at 
z ^ 6 available today. 

Aside from increasing the sample of deep spectra, there 
are several complementary observations that can improve 
these dark fraction constraints. Higher resolut i on sp ectra, 
for example the HIRES sample of I Becker et al.l (|2007l ). can 
be used to resolve flux leaks from dark pixels on smaller spa- 
tial scales. Additionally, one can use corresponding near-IR 
spectra (which we are in the process of collecting) to detect 
metal absorption corresponding to contaminating DLAs in 
the Lya/P forests, and remove their contribution from the 
dark pixel analysis. 

In addition to the observations, one could also model 
the expected contribution of the ionized IGM to the 
dark fraction. Models can be calibrated using the lat- 
est d ata: (i) analytic /numeric models of the density field 
(e.g iMiralda-Escude. Haehiielt. fc Reedbood : iTrac fc CenI 
I2OO7I : iBolton fc Becked \200^ : Cu] conservative (e.g., un- 
evolving) extrapolations of the ionization rate from lower 
redshifts, where the ionization rate is very homoge- 
neous and the Lyg forest is well resolved (e.g., ICrofg 



I2OO4I : iBolton et al.il200^ : iFaucher-Giguere et al] |2008l ): and 
(iii) conservative extrapolations of the evolution of LLSs 
I Storrie-Lombardi et al.l 19941: Stengler-Larrea etalj 19951: 



Peroux et al. L'2003'; Procha ska. O'Meara, fc Worsecl^ 



2OIOI : 



Songaila fc Co wic 2010) These approaches can be used to 



estimate the covering fraction of dark pixels expected from 
the ionized IGM and its redshift evolution. Removing this 
contamination from the total dark pixel counts leads to a 
model-dependent constraint on xm- With sufficiently high 
quality data and robust limits on the ionized IGM contribu- 
tion from models, the method presented here could be used 
to detect reionization if indeed a substantial neutral fraction 
remains at z ~ 5-6. 

Furthermore, Fig. [5] demonstrates the power of high 
dynamic range spectra from very deep observations within 
the regime where the IGM saturates absorption spectra. A 
greater dynamic range aids in detecting flux from the ionized 
IGM (except for high column density systems such as LLSs), 
but not the neutral IGM as it is far beyond the achievable 
range. Hence the rate of change of the dark fraction as a 
function of the spectral dynamic range can provide an in- 
teresting constraint on xm- Comparing this rate of change 
with the slope of an assumed density PDF could yield model- 



dependent constraints on reionization. If, for example, the 
dark fraction was observed to converge towards a constant 
value with increasing dynamic range (beyond what could be 
attributed to high column density systems), that would be 
a reionization signpost. We defer more detailed analysis to 
future work. 



8 CONCLUSIONS 

We present upper limits on the neutral hydrogen fraction at 
z ~ 5-6 derived from the simple, robust statistic of the cov- 
ering fraction of dark pixels in the Lya / /3 forests of high red- 
shift quasars. The interpretation of quasar absorption spec- 
tra is complicated by the inhomogeneous nature of reion- 
ization and the finite dynamic range of the Lyman forests: 
dark spectral regions can result either from residual H I in 
the ionized IGM, or from pre-overlap neutral patches during 
the epoch of reionization. By conservatively associating all 
dark patches with pre-overlap H I, we provide a constraint 
that is nearly model-independent. 

Using a sample of 13 z ~ 6 quasars with Keck ESI 
spectra, we constrain the neutral fraction to be xm ^ 0.2 at 
2:<5.5, rising to xm < 0.8 at z = 6.1. We find evidence of 
a break in the redshift evolution of the dark covering frac- 
tion at 2: ~ 5.5. Previous claims of a highly-ioniz ed IGM at 
these redshifts {xm ~ 10~^-10~*; e.g.. Fig. 7 in iFan et al.l 
are derived under the a priori assumption of a uni- 
form background. Given that reionization by stellar sources 
is highly inhomogeneous, this carries the implicit assump- 
tion that reionization is over, and therefore the standard 
tqp — > xhi conversion cannot be used to place direct con- 
straints on reionization, even with complete confidence in 
the density distribution of the IGM. 

At z = 6.1, more stringent constraints are provided by 
the subsample of our two deepest spectra when combined 
with conservative, albeit model-dependent, estimates of cos- 
mic variance. Specifically, we obtain Shi ^0.5. When re- 
evaluated in the context of an inhomogeneous reionization, 
the existing constraints on xui(z ~ 6) (derived from QSO 
proximity regions, damping wings in QSO and GRB spectra, 
LAE number density and clustering properties) are sensi- 
tive probes only of the early stages of reionization, when the 
Universe was mostly neutral. Additionally, they are highly 
model-dependent. Our constraint of xm{z = 6.1) < 0.5 is 
comparable to or better than these existing constraints, and 
is much more robust. 

We expect these limits to improve with a larger sample 
of deep spectra. The dark covering fraction statistic can also 
be combined with theoretical estimates of the absorption in- 
side the ionized IGM to yield a model-dependent constraint. 
Finally, the rate of change of the dark fraction as a function 
of the spectral dynamic range can provide an interesting con- 
straint on Xm- We will explore these possibilites in future 
work. 

Our results show that present-day observations of Ly- 
man forest absorption in 2 ~ 6 quasars do not rule out an 
end to reionization as late as z = 5; thus the generic state- 
ment that reionization completes by 2: ~ 6 is unjustified. In- 
deed, the tail end of reionization could stretch to 2: ~ 5 with- 
out a strong observational imprint, since the final stages of 
reionization - when the UVB became regulated by LLSs and 
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their evolution (|Furlanetto fc Mesingeij|2009l:l ICrociani et alj 
201Cl) - may have b een extended (|Furlanetto fc Qhll2005l : 



Alvarez fc Abelll2010l ). 

The goal of this work is not to promote a late reion- 
ization scenario. We simply note that it is not ruled out by 
current data, and caution against further unjustified leaps 
of interpretation. We present a direct, conservative upper 
limit on Shi that does not rely on any assumptions about 
the quasar continuum, IGM density, H II morphology or 
ionizing background fields. This can be viewed as a robust 
starting point for interpretations of high-redshift observa- 
tions and theoretical models. 
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